###################
# Post Estimation

library(rstan)
library(foreign)


# Load Datasets and combine in Data.frame
#########

  chc <- readRDS("CHES_POLICY_COMBINED_DATA_06_14_min3perc.Rds")
  nrow(unique(chc[,c("pid","year")]))
  
  c <- unique(chc$ctry_name)
  y <- c("2006","2010","2014")
  
  cd <- data.frame()
  
  for(c in unique(chc$ctry_name)){
    for(y in c("2006","2010","2014")){
      fit <- readRDS(paste("res_rozenas_",c,"_",y,".Rds",sep=""))
      d <- summary(fit,probs = c(0.01,0.025, 0.25, 0.50, 0.75, 0.975,0.99))$summary
      d <- as.data.frame(d)
      d$par <- rownames(d)
      rownames(d) <- NULL
      d$ctry_nam <- c
      d$year <- y
      cd <- rbind(cd,d)
  }}

# Check for Convergence
######
  
  # Check for Rhat
  sel_case <- unique(cd[which(cd$Rhat>=1.1),c("ctry_nam","year")]) 
  print(sel_case) 
  summary(cd$Rhat)
  
  
# Create Dataset with party estimates
#######

  s <- grep("mu\\[[0-9]+\\]|nu\\[[0-9]+\\]",cd$par)
  res_x_nu <- cd[s,]
  
  # Put PID
  res_x_nu$pid <- NA
  for(c in unique(chc$ctry_name)){
    for(y in c("2006","2010","2014")){
        dsel <- !is.na(chc[,"lr"]) 
        upids <- unique(chc$pid[dsel & chc$ctry_nam==c & chc$year==y])
        res_x_nu$pid[res_x_nu$year==y & res_x_nu$ctry_nam == c] <- rep(upids,2)
        if(length(rep(upids,2)) != length(res_x_nu$pid[res_x_nu$year==y & res_x_nu$ctry_nam == c])) cat(c,y,"\n")
    }}
  
  # Merge party names
  chall <- read.csv("1999-2014_CHES_dataset_means.csv")
  pid_nam_dta <- chall[,c("party_id","party","year")]; names(pid_nam_dta) <- c("pid","pid_nam","year")
  pid_nam_dta$pid_nam[pid_nam_dta$pid==2406] <- "TB-LNNK"
  pid_nam_dta$pid_nam <- as.character(pid_nam_dta$pid_nam)
  
  res_x_nu <- merge(res_x_nu, pid_nam_dta)
  
  
  
  # Save  
  saveRDS(res_x_nu,"res_party_rozenas.RDS")


